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D ■ Abstract 

Q 

' The measurement of the escape time of a Josephson junction might be used to detect the presence of a sinusoidal signal 
embedded in noise when standard signal processing tools can be prohibitive due to the extreme weakness of the source or to 
'yT' the huge amount of data. In this paper we show that the prescriptions for the experimental set-up and some physical behaviors 
D depend on the detection strategy. More specifically, by exploiting the sample mean of escape times to perform detection, two 
' O resonant regions are identified. At low frequencies there is a stochastic resonance/activation phenomenon, while near the plasma 
C/3 ' frequency a geometric resonance appears. Furthermore detection performance in the geometric resonance region is maximized 
^ at the prescribed value of the bias current. The naive sample mean detector is outperformed, in terms of error probability, by 
the optimal likelihood ratio test. The latter exhibits only geometric resonance, showing monotonically increasing performance 
as the bias current approaches the junction critical current. In this regime the escape times are vanishingly small and therefore 
performance are essentially limited by measurement electronics. The behavior of the likelihood ratio and sample mean detector 
for different values of incoming signal to noise ratio are discussed, and a relationship with the error probability is found. The 
f—j ' likelihood ratio test based detectors could be employed also to estimate unknown parameters in the applied input signal. As a 
O i' prototypical example we study the phase estimation problem of a sinusoidal current, that is accomplished by using the filter bank 
approach. Finally we show that for a physically feasible detector the performances are found to be very close to the Cramer- 
^ Rao theoretical bound. Applications might be found for example in some astronomical detection problems (where the all sky 
^ gravitational/radio wave search of pulsar requires the analysis of nearly sinusoidal long lived waveforms at very low signal to 
noise ratio) or to analyze weak signals in the sub-terahertz range (where the traditional electronics counterpart is difficult to 
implement). 

o 

'. I. INTRODUCTION 

^ ■ 

Threshold detection |[lll2] is based on the possibility to ascertain the presence of a signal via the transition from a metastable 
state to another Thus, the essential ingredients are: i) a two-state detector, ii) a switch from a state to the other induced by the 
signal, iii) the possibility to detect the transitions when the detector observable crosses a threshold. Under such conditions the 
original signal is transformed in a series of time-intervals, i.e. the residence times in each state. 
^ , From the physical point of view, the joint action of the applied external signal and the fluctuations induces an escape from 
• a metastable state of the detector The distribution of the escapes acquires, under very general assumptions, the shape of an 
activation law d la Arrhenius with an effective energy barrier that depends upon the fluctuations spectrum and the signal (or 
perturbation) shape H-llt]. The general idea beyond threshold detectors is that, in view of the exponential character of the 
activation law, residence times are very sensitive to small signals. Analyzing escape time series entails a loss of information 
w.rt. the direct signal observation. Conversely, it may be useful to shorten the amount of data to be analyzed. In this perspective, 
threshold detection has been examined through the lens of information theory and/or signal processing [iiSj]. 

In this paper we propose to characterize underdamped Josephson Junctions (JJ) as detectors based on the statistics of the 
escape times. It is well known that JJ are superconducting elements that can operate at extremely low temperatures (as low as 
refrigeration allows) and hence they are affected by a low intrinsic noise. Let us remark that to take advantage of the speed 
and low noise features of JJ it is necessary to analyze escape times, inasmuch the dynamics of the Josephson phase is not 
directly observable. Instead, the escape from the metastable static solution causes a voltage step, associated to the average phase 
derivative, that is actually detectable. 

Starting from the pioneering papers of Refs. i9ilTcj[]. a lot of experimental work has been performed on the role of noise in 
ac driven J J to detect the phase of the applied signal jl ij], to study quantum computation [IZl or ratchet effects ifisll . Another 
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remarkable topic is related to employ JJ as a threshold detector to characterize (up to higher order moments) weak fluctuations 
close to the quantum limit lfl4l4l7l[ (albeit there is also some controversial discussion about role of quantum noise in JJ 1, 18.1 ). In 
the ac driven quantum regime lll9ll the use of JJ is particularly appealing for signal detection since ways have been devised to 
minimize the environmental disturbances to the unavoidable quantum level |i20i] . 

JJ-based detection schemes are very interesting for the applications when the source is so weak that an intermediate amplifica- 
tion stage can introduce too much undesirable noise. Indeed SQUIDs ll2lll are sensitive to magnetic flux that is a fraction of the 
quantum flux, $o ^ 10~^^ Wb. The operating frequency range (up to the terahertz timescale) is another appealing feature of JJ. 
In this connection it can be argued that interesting applications are in the field of terahertz sensing and in the analog processing 
ll22ll . Another remarkable potential application is related to the search for gravitational wave sources where standard signal 
processing techniques (i.e. matched filter), essentially based on the intensive use of Fast Fourier Transform algorithm (FFT), are 
computationally prohibitive because applied to a large amount of data (characterized by a very small SNR) [23J. We remark that 
these examples are merely suggestive because our analysis is a proof of principle, not a carefully examination of the applications 
and of the technical problems. 

It should also be noted that a lot of theoretical problems concerning the proposed dynamical system are currently under study. 
Indeed, a resonant behavior, that can further improve detection, is expected to occur when the time scale of the external signal 
matches the fluctuation induced escape time [24.1 . Thus JJ are also a playground for topics of active research such as Stochastic 
Resonance ll25ll or Stochastic Activation ll26l - l28ll that could arise in such a scenario. 

The purpose of this paper is to provide a characterization of the JJ as a detector of a periodic signal. The aim is to show that 
nonlinear analogue processing, such as to apply the signal to a JJ, might represent a viable alternative. In fact the signal could 
be applied and processed at very high frequency (for JJ are very fast electronic devices) and with very little extra-noise from the 
nonlinear device (for JJ are superconducting elements that work also at very low temperatures). 

It has been proven that a simple statistical analysis of JJ output, such as the Sample Mean (henceforth SM) of transition times, 
could detect the presence of a sinusoid embedded in thermal noise 12911 . A more effective analysis of the transition times can be 
performed by means of appropriate methods of statistical signal processing jH^lljJ, such as the Likelihood Ratio Test (LRT). 
As expected LRT performances depend upon a suitable selection of the physical JJ parameters. The lack of analytical results 
for the escape time distribution in the so called underdamped case compels us to use extensive numerical simulations in order to 
search the JJ optimal working point in a feasible range of the relevant parameters. 

The paper is organized as follows: in Sec. II we briefly describe the physical model of the periodically driven JJ that determines 
the escape times. In Sec. Ill we describe the statistical tools used to analyze the escape times, while some technical details are 
deferred to the Appendices. In Sec. IV we present numerical results that clarify the effectiveness of the statistical analysis, both 
for signal detection and parameter estimation, namely the sinusoidal phase. Moreover we perform the JJ physical parameters 
optimization in order to further improve the detection of a perfectly known sinusoidal signal. Section V, is devoted to the 
conclusions. 



II. ESCAPE TIMES DISTRIBUTION: PHYSICAL MODEL AND MOTIVATION 

A JJ biased with a sinusoidal signal of amplitude Sq and corrupted by additive noise ^(t) is modeled by the Langevin equation 

il, 

ze ar^ HZe at 

Here C and R are the capacitance and the resistance of the JJ, respectively (we consider JJ in the underdamped regime, because 
the capacitance is not negligible). Furthermore Ic denotes the Josephson critical current, while Ib is the bias current, fc^ denote 
the Boltzmann constant and T the JJ temperature. The terms ^{t) n{t) are white Gaussian noise stochastic variables, whose 
correlators read < n{t)n{t') 25{t ~ t'), < ^{t)({t') 25{t ~ t'). 

In Eq. ([T]i thermal fluctuations \Jk\,T jR can be neglected with respect to the signal noise D 'xiT <^ DR/kb. In fact JJ can be 
cooled down at a temperature T much below the signal noise temperature, thus in the following of the paper we will assume that 
the stochastic component is dominated by the signal fluctuations of intensity D. The condition to neglect thermal fluctuations is 
favoured when the junction resistance is high, i.e. when dissipation is low, for Eq. ([T]) is based on a parallel lumped circuit model 
(see the inset of Fig. [T]i. In Sec. IV B we will show that low dissipation also favours detection, thus reinforcing the advantages 
of high R. 

Introducing the dimensionless time r — ojjt, normalized respect to the characteristic frequency (called Josephson frequency) 
ujj = [2eIc/{tiC)Y^^, and rearranging the terms, the equation in the aforementioned approximation reads 
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In Eq.© we have defined Sn = \/^jD the intensity of noise current, while the correlator in these normalized units reads 

< e(T)|(T') >= 2S{t - t') (3) 

Defining ^ — Ib/Ic the normalized bias current, a = {ujj / RIc) [fi / 2e) the normalized dissipation, e — Sq /Ic the normahzed 
signal amplitude and Em ~ {Sn /IcY 'he normalized noise intensity, Eq.(|2|i reads: 



-J^ + + sm(v3) = 7 + e sm(wT + (/jq) + V^U'^)- (4) 
A washboard potential is associated to Eq.Q that, for 7 < 1, gives rise to a barrier ll33ll 



A[/(7) = 2[Vl - 7' - 7C0S-i(7)]. (5) 

The schematic of the physics of the device is depicted in Fig. [T]as a jump over an activation barrier When the system overcomes 
the energy barrier A[/, it switches from the locked state to a running state that is associated with a finite voltage I'sT], hence it 
is possible to measure the escape time [35]. For overdamped J J the voltage step is very smooth, and it is difficult to define the 
escape from the local solution. Therefore the detection efficiency of overdamped JJ drops down. 

The main idea of signal detection is to collect the escape times to discriminate between two situations: /) the exit is caused 
by the presence of pure noise (no signal is present, Sq = 0); //) the exit is caused by the joint action of noise and a sinusoidal 
excitation (the signal is present, Sq ^ 0). Indeed escape time distributions are highly sensitive to the signal amplitude, as shown 
since the pioneering experiments 1 10]. Typical escape time Probability Density Functions (henceforth PDFs) are shown in Figs. 
|2][3]for two different sets of parameters. The random sequence of escape time is obtained by numerical integration of stochastic 
differential equation (|4]i, accumulating the first-passage-time to cross the maximum of the potential barrier (|5]l and taking into 
account the corrections of Ref. |36]. 

A solid line denotes the escape time without signal (Sq — 0) and the dashed line denotes the escape aided by a sinusoidal 
forcing (So 7^ 0). A possible signal detection strategy based on the sample mean of the escape times only measures the 
average escape time, that is essentially the slope shown in the inset of Figs. |2]|3] On the other hand. Fig. |3]makes it clear that 
the two distributions (with and without the sinusoidal signal) are, for some parameters values, very different in shape, not just 
in the mean value. This difference in shape is exploited in Section III, where it is shown that a more refined analysis can lead 
to better performances for signal detection. Moreover, in Fig. |4]it is shown the effect of the signal (initial) phase. The drastic 
change in the distribution form, while the slope remains almost constant, again demonstrates the need of a more refined analysis, 
to achieve signal phase estimation by using the escape times, as discussed in Section IV C. 

The escape time distribution not only depends upon the signal, but also on the physical JJ features. The JJ parameters that can 
be set in the experiments to achieve best detection are the electrical characteristics of the physical JJ (C, R, Ic) and the external 
bias current Ib- Moreover, if the signal is recorded, also the speed at which the recorded signal S{t) is physically applied is a 
tunable parameter A blind search among such multidimensional parameter space would be numerically prohibitive and of little 
physical significance; we therefore discuss in Section IV the physical meaning of the parameters to restrict the plausible region 
where to seek for best performances. 



III. STRATEGIES FOR JJ THRESHOLD DETECTION 

The JJ activation energy barrier defined in Eq. Q refers to an unperturbed junction, i.e. to the case 5*0 = in Eq. (|4]i. The 
addition of a deterministic sinusoidal signal of amplitude So results in an oscillation of the barrier that affects the average escape 
time. The deviations of such average escape time shown in Figs. |2][3]can be used to infer the presence of a signal. Inspection 
of Figs. |2] and |3] discloses an evident exponential decay of the escape time distribution. Such behavior has been previously 
exploited (see 112911 ) with a straightforward SM detection strategy based on the mean escape time. Indeed it can be shown that if 
the escape times are exponentially distributed in both cases, 5o = and Sq 7^ 0, then the SM detection strategy is an optimal 
strategy Bli] . In the same paper |29] the Kumar-Carrol index dKc 137] for different SNR (we recall that the SNR is related to 
the ratio Sq/Sn — Sq/ ^JZjJD) has been used as a simple (heuristic) indicator of the detector performances. We remind that 
the index dxc depends upon the operating frequency of the applied signal and the extrema are influenced by the noise intensity 
lfTTl|29ll , a peculiarity of stochastic resonance 125, 38]. 

When the PDFs are not exponential (and not explicitly known, as it is the case for underdamped JJ) we use an accurate and 
fairly common technique to determine a detection strategy and the related detector, both based on likelihood maximization. 
Such a decision criterion, based on LRT, is optimal lIsTll in a sense that is clarified below. Detectors based on LRT employ the 
knowledge of the full probability distribution of the random escape with and without the signal; therefore the shape of PDFs can 
be properly handled to improve the performances. 
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Unfortunately, the escape time distributions are not theoretically known for the system described by Eq. Even in the case 
So = the Arrhenius law is approximately valid for rare escapes |9] (in the unperturbed oscillator timescale utj), while for fast 
escapes ll39ll (that are interesting for signal analysis) only approximated analytic estimates exist. When the signal is applied the 
knowledge of the escape time distributions is even poorer, and essentially limited to the overdamped case [40]. Our solution is 
based on a semi-analytic approach, so that the analytic results are used as a guess for a numerical procedure. In the following 
we describe the proposed LRT embedded in the general framework of statistical decision theory. 



A. Statistical Decision Theory and LRT Detector 



To properly define the detection strategies, it is usual to formalize the problem as a binary hypothesis test: 

Ho : sinusoidal signal is absent 
Til : sinusoidal signal is present 

For this decision problem two different error probabilities arise: 

• the false alarm probability Pf, also called Type I error probability, i.e. the probability to decide for the hypothesis Hi 
when Ho is true; 

• the miss probability Pm, also called Type II error probability, i.e. the probability to decide for the hypothesis Ho when Hi 
is true. 

We start by considering the case in which the JJ normalized parameters (a and 7) and the normalized noise standard deviation 
Ejv are perfectly known and do not depend on the particular hypothesis in force. We also assume that the signal parameters (i.e. 
e, u! and ipo) are known under Hi hypothesis. In this setup the Neyman-Pearson lemma |30] identifies the LRT as the optimal 
detection strategy, for it minimizes, among all possible tests, the miss probability P^ at a fixed false alarm level Pj. Thus if we 
collect N escape times r = {ri},-^^^ ^j, supposed to be independent and identically distributed, the test statistic can be written 
as: 



■Hi 

h{n) > 

^_\fo{n) < 
Ho 



N 

n 



(6) 



where /o,i( ) are the PDFs of the escape times under the hypothesis Ho,i, while is a suitable threshold selected to return a 
fixed false alarm level. To simplify the computation of the statistic I©, it is possible to compare the normalized natural logarithm 
of the likelihood ratio with a threshold ( ~ \og{(')/N: 



1 ^ 



fiin. 



foin. 



Hi 

> 
< 
Ho 



(7) 



The advantage of Eq. (|7| is that the statistic A(t) can be computed as the sample mean of the random samples £ — {Ci}j^^^i ^j, 
that are obtained from the escape times via the optimal (in the Neyman-Pearson sense), non-linearity 



d = loe 



fiin) 



foin) 



(8) 



Equation (O contains the information of both PDFs /i ( • ) and /o ( • ) . Unfortunately, for underdamped JJ, an exact closed form of 
these PDFs is still unknown. We have seen in Figs. |2]and[3]that there are regimes in which both the densities follow, with a good 
approximation, an exponential law. As anticipated, in this case the SM detector is nearly optimal. Indeed it is straightforward to 
see that for exponential distributions the decision statistic in Eq. (|7| becomes: 



1 - 

Ho 



(9) 



We note that if the average escape time under hypothesis Ho is larger than the same quantity under Hi, the sign in the Eq. (|9|l 
should be reversed. 
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In other regimes we estimate both PDFs /i ( • ) and /o ( • ) using a non-parametric statistical technique such as the Kernel Density 
Estimation (KDE) 1,4 1.1 . Thus, by means of a large number of samples (5 • 10^ trials), obtained via a Monte Carlo simulation 
of the escape process of Eq. (|4|i, in both cases 5*0 = and Sq Owe retrieve a tight estimate /j (•) of the unknown PDFs /j (•). 
Further details about the KDE are given in the Appendix A. 

We are now in a position to compute the Receiver Operator Characteristic (ROC) of the test statistic, that is the plot of Pj 
vs. P„i for different values of C ll42ll - A ROC example is presented in Fig. |5] in which it is evident the trade-off between the 
two error probabilities. To simplify the performance analysis of the detector, we consider the intersection between the ROC 
and the bisector of the first quadrant angle, that is very close to the point of ROC curve with the minimum distance to the axis 
origin. In this point Pf = P„i, and we can unambiguously define the error probability P^ that is representative of the detector 
behavior The main advantage of this formulation resides in its simplicity and, in many cases, it is also a good approximation of 
the minimum Bayesian error probability when the prior probabilities of the two hypothesis Ho and "Hi are considered equal. The 
error probability Pe gives a rigorous assessment of the detector's performance, and can be related to the heuristic Kumar-Caroll 
index cIkc by an inequality, as elucidated in Appendix B. 

B. Escape Time Acquisition Strategies 

The signal to noise mixture can be applied to the JJ in different ways fl^ to acquire the escape time sequence r = {t^}^]^. 
Indeed, when the JJ switches to the running state (see Fig. ^ it is necessary to shield the junction from the signal to reset the 
static state and to apply again the signal. There is a variety of methods to reset the system. If the frequency of the signal is 
perfectly known it is in principle possible to reapply the signal always with the same initial phase. This acquisition strategy is 
called coherent detection. To apply again the signal with the same initial phase some fraction of the signal is lost waiting for the 
correct time to restart the process. In fact the correct initial times read: 

With this acquisition method the escape PDF shows a striking dependence from ipQ, as elucidated in Fig. |4] 

Another possibility, if the frequency is unknown, is to reapply the signal with any phase it might have after the reset procedure, 
and therefore with an essentially random initial value of ipo (incoherent strategy). The PDF of the escape time obtained in 
this case, essentially the average over ipo of the escape densities computed in the coherent strategy case, looses much of the 
information carried by the initial phase. 

C. Filter Bank strategy and Phase Estimation 

In the detection theory briefly summarized above, the involved parameters are supposed to be known: the LRT method can 
only be employed if the signal is supposed to have a known phase. Unfortunately in real scenarios this condition is rarely 
fulfilled and the PDF of the escape times should be considered a function of an unknown parameter vector 6. It is necessary to 
properly manage the lack of information about 8 to contrive a detection rule that minimizes the deterioration of the performances 
for unknown parameters. In this subsection we show how one can retrieve with a filter bank some information about the phase 
of the signal. The main problem arises because the signal itself might be present or not, and therefore one cannot employ a 
simple maximum likelihood procedure to determine the best guess for the parameter value. Two popular approaches to the 
detection with unknown parameter are currently used in signal processing: the Averaged LRT (ALRT) and the Generalized LRT 
(GLRT). The ALRT, based on the Bayesian Theory, consists in averaging the likelihood ratio on the unknown parameters over 
the density functions corresponding to all the values of the unknown parameter; the resulting best guess is then employed in a 
Neyman-Pearson criterion to decide the signal presence. GLRT estimates the unknown parameters via a maximum likelihood 
approach to select a single most likely phase. Such phase is substituted in the PDF expression of the LRT to decide about the 
signal presence. 

In this paper, as a paradigmatic example, we suppose that the unknown parameter is the initial phase, i.e. 0_ = ipQ. Without any 
a priori information, the ALRT approach is substantially equivalent to apply the incoherent acquisition strategy. Unfortunately 
(as is shown below) this approach looses too much information about the signal presence. Thus we focus on the GLRT, that can 
be implemented via a Filter Bank (see Fig. |6]for a pictorial scheme) . The filter bank jointly performs the estimation and the 
detection of the unknown parameter cpo- In a preliminary step we should determine a sampling of the relevant parameter space 
Wo ^Jiti such that it is suitably covered the interval [~tt,tt]. Each detector of the bank LRTi, tuned on the value i^Sq*' and 
designed as described in subsection III A, elaborates the sequence of escape times r that are acquired with a coherent strategy 
to obtain a vector of likelihood ratios {A(r, ^PQ'^)}fLi- We recall that the phase ipo, as depicted in Fig. |4] influences the PDF of 
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T only when the signal is present. Thus the maximization of the LikeUhood under Hi coincides with maximizing A(t, ipo) and 
the test (|7]i becomes 



max 

i=l,...,M 



A(r,^«)) 



Hi 

> 
< 
Ho 



(10) 



The value of (p^ ^ that maximizes A(-) constitutes an estimate of the initial phase 



ipo = argmax 

i=l,...,M 



A(r,^«)) 



(11) 



To optimize the performance of such strategy, one should design the filter bank with the appropriate choice of both the number 
M of the detector and the values ip^ K We propose to place the values ip^ ^ uniformly spaced in [— tt, tt] by sake of simplicity. 
Also, we have considered that the J J behavior for a particular phase value is not too different from another one ||29l]. The number 
M should be determined recognizing that few detectors offer a poor precision in phase estimation and in detection effectiveness, 
while too many detectors significantly increase the false alarm probability. A more precise analysis of the bank to quantitatively 
balance the two sides is out of the scope of this paper. 



IV. SIMULATION AND NUMERICAL RESULTS 



In this Section we present extensive Monte Carlo simulations that show the behavior of a detector based on a JJ. As evident in 
Eq. (|2l), LRT employs a fixed number N of escape times and it is optimal (in the Neyman-Pearson sense), under this condition. 
On the other hand the time interval Tobs to collect escape times is random and strongly depends on both the signal properties 
and the JJ parameters. Indeed, when the signal is absent (hypothesis Hq), the mean value of the escape times /io = i?[r|'Ho] , 
that is a function of 7, a and en, reads 

EiTobAno] = A^Mo- (12) 

We observe that also in presence of small signals the escape rate remains roughly unchanged fii — E[t\Hi] ~ /io- 

To perform a fair comparison under different operative conditions, we fix the mean time interval under the Hq hypothesis 
and use the corresponding number of escape times N. This strategy guarantees a mean duration time of the acquisition stage in 
any practical situation {Ho and Hi). Moreover, while the mean time of Tobs increases linearly with N, its standard deviation is 
proportional to VN, so the dispersion around the mean value becomes less significant for increasing values of N. 

Equation (fT2l l is strictly true only for the incoherent acquisition. Indeed in the coherent case we have to add to /io the mean 
time needed to assure that the incoming signal is applied to JJ with the same initial phase. So the approximate relation 

reveals that the normalized frequency lu also influences the observation time. We have proven by numerical simulations (not 
shown here) that this effect is negligible (in the interesting range of parameters), so we use the simple approximate Eq. ( fTsT i in 
the remaining part of the paper 



A. LRT Performance Improvements and Signal Frequency Analysis 

In this Subsection we investigate the dependence of the error probability Pe as a function of the normalized signal fre- 
quency cj, while in the next Subsection we carefully analyze the dependence of the performances upon the junction parame- 
ters. The distinction is somewhat artificial, inasmuch the normalized frequency uj also depends upon the junction parameters: 
uj = n/ujj = n/[2elc/ {hC)]^^^ . However, in view of some emerging physical and signal detection properties, we prefer to fo- 
cus the analysis of the signal frequency and phase in this Subsection, and postpone the analysis of the other junction parameters. 

The dependence of the detection properties upon the signal frequency is shown in Fig. |71 where we compare the SM of 
escape times (broken lines) with the more refined LRT (solid line). We also show in the inset the ROC of both the tests at a 
single frequency. In spite of the approximations used to obtain Eq. ^ the improvements of the LRT detection performances 
are significant. Along the bisector of Fig. |7]inset the probability of dismissal Pm decreases of a factor ~ lO'^ from the simple 
average to the LRT. Since the two curves are computed for the same signal duration, one can also deduce that the LRT allows 
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for a decrease of ~ 10^ of the signal length keeping the same false alarm level. It is also evident from Fig. |2]that for all driving 
frequencies the LRT overperforms the SM detector. Both strategies share a pronounced dip at the geometric resonance. Such 
frequency is not exactly ujj, for the tilted washboard potential of Eq. (|5]l exhibits a dependence of the resonant frequency upon 
the bias current of the type lf32ll : 

U:res^{l-l^f'^ (14) 

that is mirrored in the performances of the detector The resonant condition Eq. ( fT4] l is obtained through linearization of Eq. 
(|4|i for small signal amplitude e = Sq/Ic- For finite signal amplitude the oscillations explore the nonlinear part of the curvature 
that is not captured by the second order Taylor expansion behind Eq. ( fT4l i. The correction, or rectification effect, accounts for 
the deviation from Eq. ( fT4l i of the resonant frequency at finite signal amplitude, w^es — ^resi^)- In Fig- ID is shown the non 
linear J J plasma resonance curve Ii43il ; the higher order deterministic correction of the plasma frequency approaches the optimal 
working frequency for the detection. The displayed non linear plasrna^frequency is the stable branch in the regime where the 
resonance curve becomes a multivalued function of the frequency [44]. As a consequence we speculate that best detection 
performances are attained in the strongly non linear distortion regime. In our opinion this result mirrors the similar findings of 
Ref. ll45ll concerning the location of the best amplification region in the parameter space. It is noteworthy that such frequency 
optimization refinement is needed to ensure that the detector's resonance matches the signal frequency. We conclude the analysis 
of the geometric resonance affirming that both SM and LRT detectors exhibit the dependence predicted by Eq. ( fT4l l: there exists 
a suitable neighborhood of Ures that is one of the best region for detection purposes, with a small correction for the finite signal 
amplitude. 

The SM strategy reveals that a second interesting region occurs at a lower normalized signal frequency, where another reso- 
nance appears. The position of such a resonance dip depends upon the phase and the temperature |29], see Figs. |7|9f a). [T0l a). 
In fact in Ref.[|29t] it has been found that a region of optimal detection is pinpointed if the potential well barrier (tuned by 7), the 
normalized signal frequency ui (the signal frequency divided by the Josephson frequency uij) and the noise intensity D are 
connected by the relation 
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(15) 



Below the stochastic resonance frequency the PDFs (in both hypotheses) are very similar to those reported in Fig. |2l while above 
such frequency the PDF (under Hi) develops oscillations similar to those reported in Fig. |3] This explains the disappearance 
of stochastic resonance ( fTST i in the LRT detection framework, for LRT exploits the PDFs oscillations and does not deteriorate 
above the frequency (flSl l. Thus the paradoxical increase of the performances at higher noise level for the SM is solved by the 
observation that the improvement obtained at the ujsr frequency is outperformed by the choice of a more refined LRT detection 
strategy that takes into account the PDFs oscillations. In fact, see Figs. |2l the SM detector performances are always worse than 
the LRT ones, confirming the general idea that stochastic resonance is a consequence of a suboptimal detection scheme [8]. The 
practical consequence is that synergetic effects leading to stochastic resonance between noise and signal in nonlinear devices 
can only be exploited in suboptimal strategies, while in optimal detection strategies noise should only be reduced as much as the 
experimental set up allows. 



B. Physical Considerations on JJ Parameter Optimization 

In this Subsection we study the detector performances as a function of JJ parameters. Furthermore we transform the results 
obtained with simulations of the normalized Eq. (|4|i, and analyzed with the methods of Sec. Ill, into prescriptions for the 
physical parameters of an actual JJ. The system depends upon four normalized parameters that can be tuned, with an appropriate 
choice of the JJ physical parameters, to obtain the best performances. The JJ features that can be modified are: 

1. The Josephson frequency uij — [2elc/ (hC)]^^^ (and hence the capacitance C and the critical current Ic) selects the time 
scale of the device and must be chosen to maximize the performances of the device. One should compromise between two 
different requirements: fast to speed up the detection, but still slow enough to allow the available electronic to properly 
work. The critical current is typically constraint in the range IfiA < Ic < 10mA, while the capacitance is in the range 
1 — lOOOpF. The two quantities are not fully independent, for both depend upon the distance between the superconducting 
electrodes and the junction area. The critical current can also be decreased by an external magnetic field. Typically, the 
available range of uj is about 10 — lOOOGHz. 

2. The applied physical current Ib can be assumed positive (for the symmetry of the problem a negative bias value would 
just revert the phenomena) and below the critical current to have two solutions (see Eq. (|5]l), < /s < Ic- In normalized 
units the interval reads < 7 < 1. 
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3. The applied signal amplitude Sq cannot be amplified without introducing a further noise component. It is therefore 
convenient to set instead the critical current to have the more appropriate e ~ Sq/Ic normalized signal. It is worth 
noticing that the normalized intensity of the noise en = ^jD/l'^ also depends upon the critical current Ic, and can 
therefore be also tuned, while obviously the SNR cx e/y/sN = Sq / y/ujj D can not. As already mentioned the critical 
current of a Josephson junction is typically constraint in the range IfiA < 10mA, and therefore it is not possible to freely 
choose the normalized signal amplitude. Also, it is important to notice that if the critical current Ic is used to tune the 
value of the normalized signal, the system frequency ojj can still be modified via the capacitance C. 

4. The resistance R sets dissipation through the normal electron channel, parallel to the tunnel Josephson element. Dissipa- 
tion enters in Eq. (|4| through the normalized parameter a = {ujj / RIc){fi/'2e). The resistance is constrained, inasmuch 
the product IcR depends upon the material: IcR — A/2e. The energy gap A for type I traditional superconductors is in 
the range of meV. The normalized parameter a is limited in the interval < a < 1. The lower value is due to physical 
reasons (all quantities are positive), while the upper value is necessary to have a so-called hysteretic junction with two 
states (detection depends upon the possibility to detect the switch between these two states). Actual values for viable JJ 
are narrowed in the range 0.001 < a < 0.1, while the range 0.1 < a < 1 corresponds to the moderately damped regime 
where retrapping occurs, making it difficult to detect the escape from the static solution ll46ll . Finally, dissipation can be 
increased shunting the junction with an external resistor. 

One can interpret the physical parameters R, C, Ic and Is, and the corresponding normalized parameters a, u, e, and 7 in 
the following way. The normalized frequency oj = fl /ojj can be chosen to drive the JJ with the most appropriate (for signal 
detection) frequency, and this sets the ratio of the critical current Ic and capacitance C. The bias current 7 that can be varied 
through the external bias current Ib sets the potential well of the system, see Eq. (|5j, while the resistance R can be adjusted to 
have the most appropriate value of dissipation. 

As shown in Sec. IV A the normalized frequency plays a major role in the performances, see Fig. |9|a,b). A detector should 
therefore be optimized through an appropriate choice of the normalized applied frequency lu. This may be done tuning the 
Josephson angular velocity close to the signal frequency. Therefore we conclude that signals in the range of wj are suitable for 
detection with JJ, while slower or faster signals are poorly analyzed with this technique. 

Once the system is optimized in frequency, probably the most easily tunable parameter is the external bias current. Is that 
affects the normalized bias 7. In the geometric resonance neighborhood, as shown in Fig. [TOf a). there is a clear optimal point 
around 7 ~ 0.5 for the SM technique. The physical interpretation is that the average escape time is most sensitive to the external 
signal when the bias current 7, and hence the energy barrier of Eq. (|5]), is intermediate between the maximum barrier (7 = 0) and 
the minimum barrier (7 = 1). Indeed in the former case there are few escape events in the observation time, while in the latter 
case the escape is dominated by noise and the average is httle affected by the signal. It is interesting to see that for some high 
values of 7 the stochastic resonance dip around ujsb. performs better than the geometric resonance close to w^es- The LRT, as 
expected, performs better than SM for all 7 values, see Figs. [TOl' a,b). The LRT detector also shows a different behavior, namely 
the larger the bias current the more accurate the estimate, as displayed in Figs. |9lb) and [TOl' b). The qualitative explanation is 
as follows: with a more accurate analysis of the escape times it is possible to recognize the presence of the signal embedded in 
the noise even when the average is little affected by the signal itself, due to the optimal exploitation of the PDFs information - 
see Eq. (I?). It is important to notice that practical detectors cannot be realized just setting the bias current 7 = 1 for several 
reasons. First, the metastable state exists only for 7 < 1. Second, in the analysis presented here we assume that the escape time 
is measured with an infinite accuracy. In practical detectors the finite error associated to the measurements is more relevant for 
shorter escape times. We conclude that simulations suggest to use a bias current Ib as close to Ic as it is possible (i.e. 7 ~ 1), 
providing that the resulting escape times are still measurable with good accuracy. This entails that the limit of the performances 
is given by the actual electronics employed. 

The critical current might be chosen to change the normalized drive amplitude e = Sq/Ic and the noise intensity En = 
ujjD/Ic- As anticipated, the available range of the critical current is relatively narrow. However, we have found that an unifying 
parameter resumes the role of 7 (that depends upon the bias and can therefore be easily tuned) and ejv. Physical intuition 
suggests that escape time detector performances are mainly dependent on the potential barrier AU and the noise intensity en- 
A realistic guess could be that a relevant parameter is the ratio p: 



^_ A[/(7) _ 2[yr^-7Cos-i(7)] ^^^^ 

Simulations in Figs. [TT](that display the contour level of error probability as a function of AU/en and cj for two different 
values of e keeping the same SNR cx e/^/en) confirm this conjecture. Despite the nonlinear character of the system, the two 
contour plots are qualitatively similar. The main difference is due to the rectification effect (see again Fig. [8]i. We conclude 
that the analysis of parameter p leads to the same qualitative optimization recipes of parameter 7 : for the LRT strategy, the 
normalized bias current should be close to unity, while for the SM strategy optimization requires an intermediate bias value. 
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Concerning dissipation, simulations with different values a in the underdamped regime (a < 0.1) leads to the following 
result: lowering the dissipation the escape rate is increased, and therefore more events can be collected in the same time interval. 
As a consequence, the lower the dissipation the better the detector performances. So we conclude that one should try to use 
a junction resistance as large as possible to decrease dissipation. Moreover, our findings (lower dissipation favors detection) 
indicate that an external shunt to decrease the junction resistance might result in worse performances. 

Let us summarize the findings of the best JJ parameters for signal detection: 

1. The normalized frequency oj = il/ujj shows best performances around the geometrical resonance of Eq. ( fTSl ). uj ~ utres, 
for both LRT and SM. JJ are therefore best suited for signals whose frequency can be around 10 — lOOQGHz. If the 
simpler sample mean is employed, a second local optimum appears at a lower frequency, Eq. ( fTSl ), a frequency that can 
be also much smaller than the geometric resonance. 

2. The bias current for the LRT strategy should be set as close as possible to 7 = /s //c = 1 to achieve the lowest value 
of AC/, see Fig. HJb). For the SM strategy the currents should be set at an intermediate value that depends on the ratio p 
between the the energy barrier and the signal. 

3. The normalized signal eat is not an independent variable, for the relevant parameter is the ratio between the normalized 
signal and the energy barrier, Eq. ( fT6] l. This ratio should be, for the LRT strategy, as small as possible, this is automatically 
set if the bias current condition has been fulfilled . 

4. The normalized dissipation parameter a should be as low as possible. A low dissipation in fact results in faster escapes, 
thus increasing the statistics at a fixed length of the signal. The other features (in terms of optimization of the normalized 
drive frequency, bias current, and signal amplitude) are independent of dissipation. We do not show the numerical results 
for different a, for they are very similar to the displayed results obtained with a = 0.05. 

The search for optimal parameters has led to qualitative indications (e.g., the indication for high bias and matching to the 
nonlinear resonance) to select the experimental set-up without actually repeating our numerical work. The purpose of the 
previous analysis is to determine the optimal working point location in the JJ parameter space at a fixed value of the SNR. On 
the other hand, SNR strongly affects the performance of the detector evaluated through the error probability. Thus it is interesting 
to characterize the relationship between the SNR and the error probability for both SM and LRT detectors. Our ansatz is that 
the performances measured by the KC index follow the power law 



where y is a general asymptotically Gaussian decision statistics. In our setup Y is replaced by the average of the sampled 
escape times A{t) for the SM detector and by A(t) for the LRT detector, defined in Eqs. (|9|l and (O, respectively. Therefore, as 
described in Appendix B, the error probability can be expressed as 



The behavior for the optimal detection parameters has been verified via numerical simulation and the results are shown in Fig. 
fT2l^ a). It is evident that the model fit in Eq. (fTSl l is well suited for both SM and LRT detectors. More interesting is that the 
scaling law, i.e. the value of rj{Y), is different for the two detection strategies. Indeed for the SM detector we have found 
ri{A) ~ 1.63 ~ 3/2 1(29(1 . while LRT shows a scaling law ruled by ?/(A) 0.97 ^ 1, that is nearly optimal because it is 
also the behavior of the exponent for the ideal matched-filter (as can be promptly seen in Fig. (TSl b)). This striking difference 
between the SM and LRT is even more evident on the basis of the following consideration. The KC index is roughly proportional 
to y/Tobs via the coefficient B{Y), hence lowering SNR the detection time should be extended to preserve the same detector 
performance. To keep constant the quality of the detection the observation time Tots should increase as (e/ -yF/v) for the SM 

detector and as (e/y/Fiv) ^ for the LRT one. Such a different scaling gives a huge advantage to the LRT in the challenging 
case of small SNR. In Fig. \T2\ h) the comparison between the proposed LRT detection scheme and matched filter Isoil is further 
elucidated by using a suitable representation based on the KC index. In this framework, as anticipated, the 77(A) parameter 
is readily seen of the order of unity for all displayed curves. Moreover it is simple to quantify the loss in decibel (dB) of the 
detector performance, by mean of the straight lines intercepts, that for the best case is nearly 4 dB. Rephrasing this loss in terms 
of the ratio between LRT and matched filter observation time (for continuous signals) corresponds to a factor ^ 2.5. Such speed 
up factor can be achieved by frequency scaling and/or parallel VLSI implementation of the detector and JJ could turn out to be 
competitive in practical applications. In the inset of Fig. (TSlb) it is shown the dxc loss with respect to matched filter in a decibel 
scale as a function of noise parameter en- The LRT, as expected, never passes the matched filter. The minimal loss is reached 




(17) 




(18) 
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for reasonably high values of en that corresponds to the physical limit of negligible escapes towards the higher local minima. 
Finally, we have verified that inside the optimality region this function is substantially independent of the time window To^^ and 
of the ratio e / -^e/v ■ 



As anticipated in Subsection III.C, incoherent strategy may be implemented with an acquisition method that requires the 
simpler measurement setup. In Fig. [T3] three different ROC are displayed to compare the proposed detection strategies. The 
figure shows the natural behavior, i.e. that the best detection performances are obtained with the LRT detector (solid curve). A 
filter bank (GLRT coherent strategy shown in dashed curve) results in an acceptable deterioration of the detector performance. 
As anticipated, due to the huge reduction of information on the initial phase the incoherent strategy outperforms the other 
strategies. The dependence of the variance estimate (po by the GLRT detector implemented via filter bank is illustrate in Fig. [I4l 
The M — 100 sampling points are placed in ipg'' = jj{2i — 1 — M) in order to cover the interval [— tt, tt]. The sinusoidal signal 
with known parameters e = 0.1, u = 0.8 and £jv = 0.07, is applied to a JJ with a — 0.05 and 7 = 0.5. 

The unknown parameter (po = is chosen in the middle position between two sampling points, the worst case for the filter 
bank performance | |47|] . Fig. [T3]also shows (continuous curve) that the variance reduces by increasing the number N of escape 
times (the acquisition is supposed coherent). In particular the estimator variance reaches the Cramer-Rao lower bound relative 
to the escape time statistics, i.e. the asymptotic efficiency |31], for TV ^ 50 escapes. The proposed estimator is therefore able to 
effectively extract the information about the initial signal phase ipo carried by the escape times. Of course, when for very large 
values of N the expected variance becomes lower than the discretization error (27r/A/)^/12, the latter dominates and causes a 
saturation of the performances. 



We have investigated a case of detection of sinusoidal signals by means of a nonlinear device, namely an underdamped JJ. 
Some warnings are in order to clarify the limits of the approach. First, compared to traditional methods, JJ require cryogenic 
facilities to cool the superconducting electrodes. Second, the best performances are obtained at high frequency, thus requiring 
to handle microwaves without introducing distortions or undesirable noise increases. Third, we recall that a matched filter is 
the optimal method for signal detection, and therefore our proposal can be considered only when such optimal method is not 
available for technical difficulties. Forth, we underline that the improvements are obtained at the price of an heavier numerical 
work, since the escape time distribution is not theoretically known. 

We have shown that the analysis of the escape times, an experimentally accessible quantity, can be performed through an 
intuitive approach, the estimate of the average, or with a more refined likelihood-based statistic. The latter analysis leads to 
a significant improvement of the performances of the LRT with respect to SM, greatly reducing the error probability. The 
advantages of better detectors could result in solid improvements in terms of experimental time length and costs. 

We have found that the conditions that favor detection (in terms of frequency, bias current, and dissipation) depend on the 
adopted strategy (SM or LRT). Moreover, we have found that the likelihood ratio test unexpectedly leads to an optimal scaling 
law for small SNR. This result is surprising, because the signal detection has been doubly deteriorated by the insertion of a non 
linear device and the change of the observable variable (i.e. escape time instead of the original signal noise mixture). We have 
also shown that the LRT approach is scalable to detector arrays, thus allowing the signal phase estimate. However we believe 
that the most important finding is that, in some sense, the experimental setup should be tailored on the data analysis. Indeed the 
two detection approaches can not be optimized with the same normalized frequency and applied bias current. We speculate that 
the above results are fairly common to other activation detectors, because our analysis mostly relies on Arrhenius-like activation 
law, that is independent on the energy barrier details. It is therefore possible that the advantages can be exploited generalizing to 
other systems the methods we present in this work. 

Future research will be devoted to the use of signal processing tools for thermal induced and quantum assisted transition 
rates discrimination fis']. This difference is particularly difficult to detect because microwave source in J J in Macroscopic 
Quantum Tunneling regime produces multiple photon quantum transitions that can be confused with subharmonic excitations 



C. Unknown Phase Estimation and Signal Detection 



V. CONCLUSIONS 
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Appendix A: Kernel Density Estimation. 



1 ^ 



(Al) 



The Likelihood Ratio Test, presented in Eq. requires the complete knowledge of the PDFs /o( ) and /i( ) for both the 
hypothesis. As already mentioned, unfortunately there is no theoretical result that can provide this knowledge. Thus we are 
compelled to use an estimated version of these PDFs. A simple but effective strategy is to use the Kernel Density Estimation 
llirll . This technique generalizes the basic idea of histogram by using a so-called /TerneZ/Mncf/on K{-) that usually is a symmetric 
PDF. If there is a random sample X_ = {Xij-^j^ ^j, where is the sample size, the kernel estimator should be 

Nw 

i—l 

where the parameter w is the bandwidth (also called smoothing parameter). If we apply this framework to escape times, 
we immediately encounter a first difficulty. Indeed escape times are by definition positive, i.e. their PDF, under the generic 
hypothesis Hj, j G {0, 1}, has the property 

/j(i) = 0, Vt<0. (A2) 

Equation dAll l leads to an estimated PDF that does not satisfy the inequality iA2i . To deal with this issue, the following procedure 
has been applied. 

a) For a fixed hypothesis Hj, the random sample of escape times r = n] transformed via 

X = log(r). (A3) 
Thus we deal with the random sample X = {-'^iljgji ^v] ^^'^^ '-^'^ assume every value on the real axis. 

b) By means of Eq. ( lAll ). an estimated PDF gj (x) is computed. 



c) The PDF fj (t) is obtained from gj (x) via 



i- (t) = ^«^, t>0. (A4) 



The procedure described above has to be applied twice, i.e. under both the hypothesis Hj, to obtain the estimates /o(') ^nd /i(-) 
PDF. Then, when we have to decide about the presence of a signal via an independently generated random sample of escape 
times, we can compute the statistic in Eq. dTji as: 



1 ^ 

1=1 



5i(log(T,)) 



?o(log(T,)) 



•Hi 

> 



< C, (A5) 



Ho 

The procedure is completely specified if in Eq. dAlb one chooses both the kernel function and the bandwidth. We have used the 
kernel defined by a standard Gaussian density, i.e. 

This kernel is rewarding, for the estimated density is smooth and there are no subset on the real axis with zero density. The latter 
point is essential in the likelihood ratios to avoid singularities. The drawback is in a little loss in estimation performances, that 
is negligible for the used sample size (^ 5 • 10^). The smoothing parameter selection, instead is an hard task. Indeed it is known 
(see MM) that the optimal choice has the form 

w = C{K,f)N~^/^ (A6) 

where C{K, f) is a constant that depends on the used kernel K{-) and, unfortunately, on the same density /(•) to be estimated. 
Thus we have used a simple closed form obtained when both K{-) and /(•) are Gaussian that is 

w=i^] dN-^'\ (A7) 

where a is the sample standard deviation of the dataset. This choice is near-optimal for unimodal densities, that is the case 
when the Ho hypothesis is in force, while it leads to some over-smooth for multimodal densities that arise under "Hi hypothesis. 
Also in this case the large sample used helps us make the over smoothing negligible. Moreover, the slight over-smoothing 
introduced under Hi can only worsen the LRT test performance, for the PDF oscillations, that contain the largest part of 
additional information with respect to the sample mean, are underestimated. 
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Appendix B: Kumar-Carrol index. 



In this Appendix we introduce the general concept of KC index and its connection with a more effective parameter: the error 
probability. In the contest of decision theory , the main purpose is to discriminate between two hypothesis: signal is present (Hi) 
vs signal is absent CHo). Let Y the selected decision statistics, expressed as function of the escape times vector r, that has to be 
compared with a suitable threshold Moreover, suppose that and fJ,o{Y) are the averages of the decision statistic under 

the two hypothesis (subscripts 1 and refer to the presence and absence of the signal, respectively) while the corresponding 
standard deviations are denoted with ai{Y),(Jo{Y). The KC index can be accordingly defined as 

. I I . (Bi, 
W(y) + »g(y)) 

In this paper both the decision statistics can be modelled as the sample mean of a suitable random variable, i.e. r for SM detector 
and C for LRT one (see Eqs. Q and ^ respectively). In this case both the statistics, for large sample size, are asymptotically 
normal due to the Central Limit Theorem |49]. Therefore the error probability (defined above as the value of the ROC 
Pjn = Pf) can be expressed as 




f~A(YWdKc(Y)\ 

— ^ 



(B2) 



where A{Y) — 2 | (Jq{Y) — o-i{Y) \ / \ o-q{Y) + (Ji{Y) \. By inspection it can be shown that Eq. ( IB2l i is a decreasing function 
of A, and therefore neglecting the difference among standard deviations (if it exists) it is possible to retrieve an upper bound of 
Pe that is only function of dxc, i-e- 

The inequality ( IB3b clarifies the heuristic character of the KC index as an indicator of the detector performance. The concept 
above introduced can be applied to SM and LRT detectors, by setting Y := A{z) and Y := A(r), respectively. Under the 
hypothesis that a power law relationship between KC index and SNR exists, as exemplified by 

dKc{Y) ^ A{Y) i^-^j . (B4) 

The error probability can be expressed as 

where B{Y) = \J ^ + -^^j^ / (2a/2). This latter equation has been used to interpret the relationship between Pg and 
e/ ^Jen in Sec. IV, see Fig. [12] 
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FIG. 2: Distribution of the escape times witli thie sinusoidal signal (dashed line) and without (continuous line). The main effect of the signal 
for this range of the parameters is to change the slope of the distribution. The inset shows the asymptotic fit with exponential distributions. 
The resulting decay rate in presence of (without) the signal is 0.069 (resp. 0.055). Parameters of the simulations are: 7 — 0.5, a = 0.05, 
£jv = 0.07. Moreover, when the signal is present, e = 0.1, ipo — and oj = 0.035. 
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FIG. 3: Distributions of the escape times witli tlie sinusoidal signal (dashed line) and without (continuous line). The main effect of the signal 
for this range of the parameters is to introduce oscillations on the escape time distribution. The inset shows the asymptotic fit with exponential 
distributions. The resulting decay rate in presence of (without) the signal is 0.066 (resp. 0.055). Parameters of the simulations are: 7 = 0.5, 
a = 0.05, ejv = 0.07. Moreover, when the signal is present, e = 0.1, c^o = and w = 0.8. 




FIG. 4: Distributions of the escape times with a signal with two initial phases, ipo = and (fio — tt. The inset shows the asymptotic fit with 
exponential distributions, it is evident that the exponential decay rate is almost the same. Indeed the resulting decay rate for = (i^o = w) 
is 0.067 (0.068). The other parameters of the simulations are: 7 = 0.5, a = 0.05, en = 0.07, e = 0.1 and w = 0.8. 



16 




FIG. 5: Typical ROCs of a JJ based detector under the hypothesis of complete parameter knowledge for different values of dxc (i.e. of e). 
Other relevant simulation parameters are: 7 = 0.3, a = 0.05, en = 0.07. Moreover, when the signal is present, (^0 = and w = 0.8. 




FIG. 6: Schematic of the filter bank. The JJ block indicates the physical device giving the escape time, while the Di blocks represent the fiher 
bank in which the elements are tuned on different values of the unknown parameter (in this case the initial phase ipo). The maximization 
block (MAX) output is compared with a suitable threshold C to perform the decision test while arg max gives the maximum likeUhood phase 
estimation. 
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FIG. 7: Detection probability as a function of the applied signal frequency. Parameters of the system are 7 — 0.8, a = 0.05, en ~ 0.0175, 
£ — 0.05 and ifo — 0. The simulations are performed setting the mean observation time under Ho, E[Tobs\'Ho] = 2000. In the inset the 
ROCs computed for both the detectors at a frequency close to geometric resonance a;,.es ~ cj = 0.7. 




FIG. 8: Plot of frequency u where a minimum of Pe. occurs as a function of the bias current 7 for LRT detection strategy. Dashed line refer 
to a signal with e = 0.1 and a noise with ejv = 0.07, while dotted line to e = 0.05 and a noise with en ~ 0.0175 (see also Fig. llOt . For 
sake of comparison we also show (continuous line) the resonant frequency of Eq. il4l corresponding to vanishingly small e and the non linear 
deterministic resonance curve (dots) in the case e — 0.05. Other relevant parameters of the system are; a = 0.05 and (po = 0. The simulations 
are performed setting the mean observation time under Ho, E[Tobs\'Ho] — 2000. 
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FIG. 9: Multiple plot of the error probability Pe as a function of the applied signal frequency cj and the bias current 7. (a): SM detector, (b) 
LRT detector. Along the arrow bias 7 decreases from 0.9 to 0.4 with step —0.1. Other relevant parameters of the system are: a = 0.05, 
sn = 0.0175, e = 0.05 and tpo = 0. The simulations are performed setting the mean observation time under Ho, ElTobslHo] = 2000. 
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FIG. 10: Contour plot of the error probability Pe as a function of the applied signal frequency u) and the bias current 7. (a): SM detector, (b) 
LRT detector. Thick line indicates in (b) the location of minimum Pe probabiUty. Other relevant parameters of the system are: a = 0.05, 
ejv = 0.0175, e = 0.05 and (po = 0. The simulations are performed setting the mean observation time under Ho, ElTotslHo] = 2000. 
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FIG. 1 1: Contour plot of the error probability Pe as a function of the applied signal frequency uj and the noise-normalized barrier AU/en- (a) 
e = 0.1 and en = 0.07, (b) s = 0.05 and sn ~ 0.0175. Other relevant parameters of the system are: a = 0.05 and (po = 0. The simulations 
are performed setting the mean observation time under Ho, E\Toba\Ho\ = 2000. 
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FIG. 12: In the part (a) we display the error probabihty Pe as a function of e/^/en (related to SNR). Disks are obtained by numerical 
simulations, continuous lines are computed by fitting the data using the KC approximation in Eq. l lISI l. The parameters used in the simulations 
are: 7 = 0.3, a — 0.05, £jv = 0.07. Moreover, when the signal is present, tpo = and uj — 0.8, close to the optimal detection point ujres- 
Below in (b) the KC index (reported in dB) of proposed LRT strategy (for different parameters 7 = 0.3, 0.9 and cj near to the optimal detection 
frequency) vs. the same SNR related parameter is compared with the KC index of the Matched Filter (best available strategy). In the inset the 
dKC loss w.r.t. matched filter (in a dB scale) for 7 = 0.9 is displayed as a function of ejv (e/^ejv is fixed to —20 dB). In all cases simulations 
are performed with ElTotslHo] = 2000. 
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FIG. 13: ROCs computed for three different detectors: LRT detector based on incoherent acquisition of escape times (dot dashed line); filter 
bank coherent acquisition (GLRT) detector (dashed line) and LRT detector based on coherent data (continuous line). For all strategies we 
ensure the same mean observation time under Ho (i.e. E[Tobs \Ho] w 1500). Other relevant parameters are: 7 = 0.5, a = 0.05, £n = 0.07. 
Moreover, when the signal is present, e = 0.1, i^o = and oj = 0.8, close to the optimal detection point ojres- 




FIG. 14: Variance of the maximum likelihood estimate of (po as a function of the number of escape times. Continuous line is obtained by 
Monte Carlo method, dashed line is the theoretical Cramer-Rao lower bound. In the inset a typical probability mass function estimate (with 
10^ trials) of the unknown parameter tpo by using only N = 10 samples. The signal is generated with an initial phase ipo = 0. Other relevant 
parameters are: 7 = 0.5, a = 0.05, en = 0.07, e = 0.1 and to = 0.8, close to the optimal detection point Wres- 



